# =============================================================================================
# - THIS FILE REPLICATES TABLES 1--5B and FIGURE 2 OF THE PAPER
# =============================================================================================

# =============================================================================================
# LOAD DATA
# =============================================================================================

# load packages
library(data.table); library(magrittr); library(fastmatch); library(pbapply); library(xts)
library(timeDate); library(bizdays); library(fasttime)
library(xtable); library(stargazer); library(haven); library(lfe)
library(ggplot2); library(cowplot); library(grid); library(ggridges)
library(bife); library(scales);library(texreg); library(alpaca)

# color-blind palette
cbp <- c("#000000", "#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2", "#D55E00", "#CC79A7",
         "#888888")
bdf_pale <- rgb(0.192, 0.549, 0.906, alpha = 0.35)
bdf <- rgb(0.192, 0.549, 0.906, alpha = 1)

# set locale
Sys.setlocale(locale = 'US')

# load user-defined function
source('desc.R')

load_rmetrics_calendars(seq(2002, 2023, 1))

# load
load(file = 'replica.RData')

# =============================================================================================
# TABLE 1
# =============================================================================================

# summary stats
sumtab = as.data.frame(t(sub_event_df[,
                                      sapply(.SD, function(x) { as.numeric(descsubfull(x)) }),,
                                      .SDcols = c('mat_pre',
                                                  'size_hist_pre',
                                                  'age_pre',
                                                  'coupon_pre',
                                                  'rating_num_pre',
                                                  'duration_pre',
                                                  'ys_pre',
                                                  'dts',
                                                  'calliv_mkt_pre',
                                                  'event_call',
                                                  'mat_diff_mkt',
                                                  'size_diff_mkt',
                                                  'age_diff_mkt',
                                                  'coupon_diff_mkt',
                                                  'duration_diff_mkt',
                                                  'ytm_diff_mkt',
                                                  'mve_bln_pre',
                                                  'bl1_pre',
                                                  'pm_pre',
                                                  'cash_pre',
                                                  'ret_equity_prior')]))

colnames(sumtab) = names(descsubfull(sub_event_df$ret))
rownames(sumtab) = c('Maturity, years',
                     'Size, \\$ mn',
                     'Age, years',
                     'Coupon rate, \\%',
                     'Rating',
                     'Duration, years',
                     'Yield spread, \\%',
                     'Return volatility, \\%',
                     'Money-ness, \\%',
                     'Called (dummy)',
                     '$\\Delta$Maturity, years',
                     '$\\Delta$Size, \\$ mn',
                     '$\\Delta$Age, years',
                     '$\\Delta$Coupon rate, \\%',
                     '$\\Delta$Duration, years',
                     '$\\Delta$Yield spread, \\%',
                     'Market cap., \\$ bn',
                     'Book leverage, \\%',
                     'Profit margin, \\%',
                     'Cash ratio, \\%',
                     'Equity return, \\%')

digimat <- cbind(rbind(matrix(nrow = nrow(sumtab), ncol = 8, 2)),
                 matrix(nrow = nrow(sumtab), ncol = 1, 0))
digimat[c(12), 1:(ncol(digimat)-1)] <- 0

addtorow = list()
addtorow$pos = list(0, 10, 16)
addtorow$command = c("\\hline \\multicolumn{9}{c}{(A) Callable bonds characteristics}  \\\\\n",
                     "\\hline \\multicolumn{9}{c}{(B) Difference to duration- and rating-matched portfoio}  \\\\\n",
                     "\\hline \\multicolumn{9}{c}{(C) Issuer characteristics}  \\\\\n")

print(xtable(sumtab, auto = F,
             align = c('l|', rep('Y', '8')),
             digits = digimat),
      file = 'T1.tex', floating = F,
      include.colnames = T,
      tabular.environment = 'tabularx',
      width = '\\textwidth', 
      include.rownames = T, size = 'footnotesize', sanitize.text.function=function(x){x},
      hline.after = c(nrow(sumtab)), add.to.row = addtorow)

# =============================================================================================
# TABLE 2A
# =============================================================================================

# variable names
varnms <- c('$R$', '$XR$', '$XRS^{b}$', '$XRS^{n}$')
retvec <- c('ret', 'xsret', 'xsret_exmkt', 'xsret_exnclbl')
daretvec <- paste0('da', retvec)

m1 = lapply(as.list(c(retvec[-c(3)], daretvec[-c(3)])), function(x) {
  
  lm(formula(paste0(x, '~ 1')), sub_event_df)
  
})

m2 = lapply(as.list(c(retvec[-c(3)], daretvec[-c(3)])), function(x) {
  
  lm(formula(paste0(x, '~ as.factor(per)-1')), sub_event_df)
  
})

mod1 <- capture.output(
  stargazer(out = 'T2a.tex', m1,
            align = T,
            font.size = "footnotesize", no.space = T, omit.stat = c('ser', 'f', 'adj.rsq', 'rsq', 'n'),
            covariate.labels = c('Full sample'),
            float = F, intercept.top = T, intercept.bottom = F, digits = 2,
            dep.var.labels.include = F, model.numbers = F, model.names = F,
            dep.var.caption = "",
            #column.labels = c("All", yrs),
            column.labels = c('$~~R~~$', '$~XR~$', '$XRS^{n}$',
                              '$~~r~~$', '$~xr~$', '$xrs^{n}$'),
            notes.align = 'l', column.sep.width = "0pt")
)
#mod1 <- sub("lcccccc", "lccc|ccc", mod1)
mod1[10] <- paste(mod1[10], "\\hline")

mod2 <- capture.output(
  stargazer(out = 'T2a.tex', m2,
            align = T,
            font.size = "footnotesize", no.space = T, omit.stat = c('ser', 'f', 'adj.rsq', 'rsq', 'n'),
            covariate.labels = c('Pre-GFC', 'GFC (2008/09)', 'Post-GFC'),
            float = F, intercept.top = T, intercept.bottom = F, digits = 2,
            dep.var.labels.include = F, model.numbers = F, model.names = F,
            dep.var.caption = "",
            #column.labels = c("All", yrs),
            column.labels = c('$~~R~~$', '$~XR~$', '$XRS^{n}',
                              '$~~r~~$', '$~xr~$', '$xrs^{n}'),
            notes.align = 'l', column.sep.width = "0pt")
)
#mod2 <- sub("lcccccc", "lccc|ccc", mod2)

mod <- c(mod1[1:13],
         mod2[11:16],
         mod2[18:22])

mod[7] <- "\\begin{tabular}{@{\\extracolsep{0pt}}lD{.}{.}{4} D{.}{.}{4} D{.}{.}{4} | D{.}{.}{4} D{.}{.}{4} D{.}{.}{4} } "
mod[10] <- " & \\multicolumn{1}{c}{$~~R~~$} & \\multicolumn{1}{c}{$~XR~$} & \\multicolumn{1}{c|}{$XRS^{n}$} & \\multicolumn{1}{c}{$~~r~~$} & \\multicolumn{1}{c}{$~xr~$} & \\multicolumn{1}{c}{$xrs^{n}$} \\\\  \\hline"

writeLines(mod, 'T2a.tex')

# =============================================================================================
# TABLE 2B
# =============================================================================================

m <- vector('list', 12)

m[[1]] <- felm(ret ~ cm + mat + size_hist | rating + month | 0 | issuer, mtrace)
m[[2]] <- felm(ret ~ per:cm + mat + size_hist | rating + month | 0 | issuer, mtrace)
m[[3]] <- felm(ret ~ cm + mat + size_hist | rating + month | 0 | issuer, mtrace[rating_num <= 10])
m[[4]] <- felm(ret ~ per:cm + mat + size_hist | rating + month | 0 | issuer, mtrace[rating_num <= 10])
m[[5]] <- felm(ret ~ cm + mat + size_hist | rating + month | 0 | issuer, mtrace[rating_num > 10])
m[[6]] <- felm(ret ~ per:cm + mat + size_hist | rating + month | 0 | issuer, mtrace[rating_num > 10])

m[[7]] <- felm(xsret_alt ~ cm + mat + size_hist | rating + month | 0 | issuer, mtrace)
m[[8]] <- felm(xsret_alt ~ per:cm + mat + size_hist | rating + month | 0 | issuer, mtrace)
m[[9]] <- felm(xsret_alt ~ cm + mat + size_hist | rating + month | 0 | issuer, mtrace[rating_num <= 10])
m[[10]] <- felm(xsret_alt ~ per:cm + mat + size_hist | rating + month | 0 | issuer, mtrace[rating_num <= 10])
m[[11]] <- felm(xsret_alt ~ cm + mat + size_hist | rating + month | 0 | issuer, mtrace[rating_num > 10])
m[[12]] <- felm(xsret_alt ~ per:cm + mat + size_hist | rating + month | 0 | issuer, mtrace[rating_num > 10])

mod1 <- capture.output(
  stargazer(out = 'T2b.tex', m[c(1,3,5,7,9,11)],
            align = T,
            font.size = "footnotesize", no.space = T, omit.stat = c('ser', 'f', 'adj.rsq', 'rsq', 'n'),
            covariate.labels = c('Full sample'),
            float = F, intercept.top = T, intercept.bottom = F, digits = 2,
            dep.var.labels.include = F, model.numbers = F, model.names = F,
            dep.var.caption = "",
            column.separate = c(3,3),
            column.labels = c("$\\tilde{R}$", '$\\tilde{XR}$'),
            omit = c('mat', 'size_hist'),
            order = c('^perPre', '^perGFC', '^perPost'),
            notes.align = 'l', column.sep.width = "0pt")
)

mod2 <- capture.output(
  stargazer(out = 'T2b.tex', m[c(2,4,6,8,10,12)],
            align = T,
            font.size = "footnotesize", no.space = T, omit.stat = c('ser', 'f', 'adj.rsq', 'rsq', 'n'),
            covariate.labels = c('Pre-GFC', 'GFC', 'Post-GFC'),
            float = F, intercept.top = T, intercept.bottom = F, digits = 2,
            dep.var.labels.include = F, model.numbers = F, model.names = F,
            dep.var.caption = "",
            omit = c('mat', 'size_hist'),
            order = c('^perPre', '^perGFC', '^perPost'),
            column.separate = c(3,3),
            column.labels = c("$\\tilde{R}$", '$\\tilde{XR}$'),
            notes.align = 'l', column.sep.width = "0pt")
)

mod <- c(mod1[1:10],
         ' & \\multicolumn{1}{c}{$\\text{All}$} & \\multicolumn{1}{c}{IG} & \\multicolumn{1}{c|}{HY} & \\multicolumn{1}{c}{All} & \\multicolumn{1}{c}{IG} & \\multicolumn{1}{c}{HY} \\\\',
         ' \\hline \\\\[-1.8ex]',
         mod1[11:13],
         mod2[11:16],
         mod2[18:22])

#mod <- sub("lcccccc", "lccc|ccc", mod)
mod[7] <- "\\begin{tabular}{@{\\extracolsep{0pt}}lD{.}{.}{4} D{.}{.}{4} D{.}{.}{4} | D{.}{.}{4} D{.}{.}{4} D{.}{.}{4} } "

writeLines(mod, 'T2b.tex')

# =============================================================================================
# TABLE 5A and 5B
# =============================================================================================

# returns and volumes by money-ness
retvec_ext <- c(c(retvec, daretvec), 'ret10k', 'xsret10k', 'dnetvlm', 'netvlm_pre')
retvec_ext_nms <- c(gsub("\\$", "", c(varnms, tolower(varnms))),
                    "R^{ex10}", "XR^{ex10}",
                    "\\Delta ADNV^{\\text{post-pre}}", "ADNV^{\\text{pre}}")

# make a joint table with common effects for all bonds
sub_event_df[, itm_mkt_pre_reg := factor(itm_mkt_pre, levels = c('ITM', 'ATM', 'OTM'))]
sub_event_df[, itm_mkt_adj_pre_reg := factor(itm_mkt_adj_pre, levels = c('ITM', 'ATM', 'OTM'))]


  for (j in 1:2) {
    
    if (j == 1) {
      
      df <- sub_event_df
      nmf <- 'a'
      
    } else {
      
      df <- sub_event_df[per=='postgfc']
      nmf <- 'b'
      
    }
    
    vrbl <- retvec_ext[2]
    vrbl_nm <- retvec_ext_nms[2]
    m <- vector('list', length = 1)
    
    keypart <- '~ itm_mkt_pre_reg + surprise_dummy:itm_mkt_pre_reg'
    
    if (grepl('exmkt', vrbl)) {
      
      mod <- '+ size_diff_mkt + mat_diff_mkt + coupon_diff_mkt'
      toomit <- c('size_diff_mkt', 'mat_diff_mkt', 'coupon_diff_mkt')
      
    } else if (grepl('nclbl', vrbl)) {
      
      mod <- '+ size_diff + mat_diff + coupon_diff'
      toomit <- c('size_diff', 'mat_diff', 'coupon_diff')
      
    } else {
      
      mod <- '+ size_hist_pre + mat_pre + coupon_pre'
      toomit <- c('size_hist_pre', 'mat_pre', 'coupon_pre')
      
    }
    
    m[[1]] <- felm(formula(paste(vrbl, keypart,
                                 '- 1|0|0|issuer')), cmethod = 'reghdfe', 
                   data = df)
    
    m[[2]] <- felm(formula(paste(vrbl, keypart, mod,
                                 '|rating_num_pre|0|issuer')), cmethod = 'reghdfe', 
                   data = df)
    
    m[[3]] <- felm(formula(paste(vrbl, keypart, mod,
                                 '|year+issuer+rating_num_pre|0|issuer')), cmethod = 'reghdfe', 
                   data = df)
    
    m[[4]] <- felm(formula(paste(vrbl, keypart, mod,
                                 '|issuerxyear+ rating_num_pre|0|issuerxyear')), cmethod = 'reghdfe', 
                   data = df[issuerxyear %in% largeissuers])
    
    m[[5]] <- felm(formula(paste(vrbl, keypart, mod,
                                 '+dts|issuerxyear+ rating_num_pre|0|issuerxyear')), cmethod = 'reghdfe', 
                   data = df[issuerxyear %in% largeissuers])
    
    m[[6]] <- felm(formula(paste(vrbl, keypart, mod,
                                 '+dts+ret_equity|issuerxyear+ rating_num_pre|0|issuerxyear')), cmethod = 'reghdfe', 
                   data = df[issuerxyear %in% largeissuers])
    
    m[[7]] <- felm(formula(paste(vrbl, keypart, mod,
                                 '+dts+ret_equity+dnetvlm|issuerxyear+ rating_num_pre|0|issuerxyear')), cmethod = 'reghdfe', 
                   data = df[issuerxyear %in% largeissuers])
    
    mod <- capture.output(
      
      stargazer(out = paste0('T5', nmf, '.tex'),
                align = T,
                m,
                font.size = "footnotesize", no.space = T, omit.stat = c('ser', 'rsq'),
                covariate.labels = c('$\\mathbb{1}^{\\text{Surprise}}\\times\\mathbb{1}^{\\text{ITM}}$',
                                     '$\\mathbb{1}^{\\text{Surprise}}\\times\\mathbb{1}^{\\text{ATM}}$',
                                     '$\\mathbb{1}^{\\text{Surprise}}\\times\\mathbb{1}^{\\text{OTM}}$',
                                     '$\\text{Return volatility}$',
                                     'Equity return',
                                     '$\\Delta ADNV$',
                                     '$\\mathbb{1}^{\\text{ITM}}$',
                                     '$\\mathbb{1}^{\\text{ATM}}$',
                                     '$\\mathbb{1}^{\\text{OTM}}$'),
                float = F, intercept.top = T, intercept.bottom = F, digits = 2,
                omit = toomit,
                model.numbers = T, model.names = F,
                dep.var.caption = paste0("$\\text{Dependent variable:~}", vrbl_nm, '$'),
                dep.var.labels.include = F,
                order = c(':',
                          'dts',
                          'ret_equity',
                          'dnetvlm',
                          '^itm_mkt_pre_regITM$',
                          '^itm_mkt_pre_regATM$',
                          '^itm_mkt_pre_regOTM$'),
                add.lines = list(c("Bond controls", rep('\\multicolumn{1}{c}{$\\text{NO}$}', 1), rep('\\multicolumn{1}{c}{$\\text{YES}$}',6)),
                                 c("Year, Issuer FEs", rep('\\multicolumn{1}{c}{$\\text{NO}$}', 2), rep('\\multicolumn{1}{c}{$\\text{YES}$}',1), rep('\\multicolumn{1}{c}{$\\text{NO}$}',4)),
                                 c("Issuer$\\times$Year FE", rep('\\multicolumn{1}{c}{$\\text{NO}$}', 3), rep('\\multicolumn{1}{c}{$\\text{YES}$}',4))),
                notes.align = 'l', column.sep.width = "0pt")
      
    )
    
    mod[7] <- "\\begin{tabular}{@{\\extracolsep{0pt}}lD{.}{.}{4} D{.}{.}{4} D{.}{.}{4} D{.}{.}{4} D{.}{.}{4} D{.}{.}{4} D{.}{.}{4} } "
    
    writeLines(mod, paste0('T5', nmf, '.tex'))
    
  }

# =============================================================================================
# TABLE 4
# =============================================================================================

m = vector('list')
sub_event_df[, called := as.numeric(event_call)*100]
sub_event_df[, calliv_mkt_pre_itm := calliv_mkt_pre*as.numeric(itm_mkt_pre%in%c('ATM', 'ITM'))]

extradf <- sub_event_df[, .(issuer, cusip_id, year, issuerxyear, event_date, called, xsret, event_id,
                            event_call, calliv_mkt_pre, busy, dts, ysvol_pre, ytmvol_pre, logpretvol_pre,
                            size_hist_pre, rating_num_pre, mat_pre, rating, ret_equity, ret_equity_prior,
                            cash_pre, ysdelta_pre,
                            mve_bln_pre, bl1_pre, bl2_pre, ml_pre, atan_pre, pm_pre, tpr_pre, capexr_pre)]
colstosd <- c('calliv_mkt_pre', 'dts', 'ysvol_pre', 'ytmvol_pre', 'logpretvol_pre',
              'size_hist_pre', 'mat_pre', 'ret_equity', 'ret_equity_prior', 'cash_pre', 'ysdelta_pre',
              'mve_bln_pre', 'bl1_pre', 'bl2_pre', 'ml_pre', 'atan_pre', 'pm_pre', 'tpr_pre', 'capexr_pre')

extradf[, (colstosd) := lapply(.SD, function(x){x/sd(x, na.rm = T)}), .SDcols = colstosd]
extradf[, HY := 0]
extradf[rating_num_pre > 10, HY := 1]

setkey(extradf, issuer, event_date)
extradf[, nbevent := seq(1, .N, 1), .(issuer, cusip_id)]
extradf[, (paste0('lag_', 1:2, '_IV')) := shift(calliv_mkt_pre, n = 1:2, type = 'lag'), .(cusip_id)]
extradf[!is.na(lag_1_IV), ivup := 0]
extradf[calliv_mkt_pre - lag_1_IV > 0, ivup := 1]

m = vector('list')

m[[1]] <- glm(event_call ~ calliv_mkt_pre,
              data=extradf, family = binomial('logit'))

m[[2]] <- glm(event_call ~ calliv_mkt_pre + dts,
              data=extradf, family = binomial('logit'))

m[[3]] <- feglm(event_call ~ calliv_mkt_pre + dts|year|issuer,
                data=extradf, family = binomial('logit'), control = list(drop.pc = F))

m[[4]] <- feglm(event_call ~ calliv_mkt_pre + dts|year+rating_num_pre|issuer,
                data=extradf, family = binomial('logit'), control = list(drop.pc = F))

m[[5]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre|year+rating_num_pre|issuer,
                data=extradf, family = binomial('logit'), control = list(drop.pc = F))

m[[6]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre+
                  mve_bln_pre + bl1_pre |year+rating_num_pre|issuer,
                data=extradf, family = binomial('logit'), control = list(drop.pc = F))

m[[7]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre+
                  mve_bln_pre + ml_pre |year+rating_num_pre|issuer,
                data=extradf, family = binomial('logit'), control = list(drop.pc = F))

m[[8]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre+
                  mve_bln_pre + bl1_pre + pm_pre + tpr_pre + cash_pre|year+rating_num_pre|issuer,
                data=extradf, family = binomial('logit'), control = list(drop.pc = F))

m[[9]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre+
                  mve_bln_pre + bl1_pre + pm_pre + tpr_pre + cash_pre|year+rating_num_pre+issuer|issuer,
                data=extradf, family = binomial('logit'), control = list(drop.pc = F))

m[[10]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre+
                   mve_bln_pre + bl1_pre + pm_pre + tpr_pre+ cash_pre +
                   busy|year+rating_num_pre+issuer|issuer,
                 data=extradf, family = binomial('logit'), control = list(drop.pc = F))

m[[11]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre+
                   mve_bln_pre + bl1_pre + pm_pre + tpr_pre+ cash_pre+
                   busy+ ivup|year+rating_num_pre|issuer,
                 data=extradf, family = binomial('logit'),
                 control = list(drop.pc = F))

m[[12]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre+busy|issuerxyear|issuerxyear,
                 data=extradf[issuerxyear %in% largeissuers], family = binomial('logit'),
                 control = list(drop.pc = F))

m[[13]] <- feglm(event_call ~ calliv_mkt_pre + dts+size_hist_pre+mat_pre+busy+ivup|issuerxyear|issuerxyear,
                 data=extradf[issuerxyear %in% largeissuers], family = binomial('logit'),
                 control = list(drop.pc = F))

# extract ribust SEs and PVs
sel <- list('SE' = vector(mode = 'list', length = length(m)),
            'PV' = vector(mode = 'list', length = length(m)))

for (i in 1:length(m)) {
  
  if (i %in% c(3:11)) {
    
    A <- summary(m[[i]], type = 'clustered', cluster = ~ issuer)$cm
    
  } else if (i %in% c(12,13)) {
    
    A <- summary(m[[i]], type = 'clustered', cluster = ~ issuerxyear)$cm
    
  } else {
    
    A <- summary(m[[i]])$coefficients
    
  }
  
  sel$SE[[i]] <- as.numeric(A[,2])
  sel$PV[[i]] <- as.numeric(A[,4])
  
}

texreg(file = 'T4.tex',
       l = m,
       custom.header = list("Dependent variable: $\\mathbf{P}_i\\left[\\text{Called}\\right]$" = 1:length(m)),
       custom.model.names = paste0('(',1:length(m),")"),
       custom.coef.names = c('Intercept',
                             'Intrinsic value (IV)',
                             'Return volatility',
                             'Outstanding amount', 'Maturity',
                             'Equity value', 'Book leverage', 'Market leverage',
                             'Profit margin', 'Total payout ratio', 'Cash ratio',
                             'Fiscal year-end (dummy)',
                             'IV > IV @ prev. notice (dummy)'),
       #reorder.coef = c(1,2,4,3,5:15),
       override.se = sel$SE,
       override.pvalues = sel$PV,
       include.deviance = F,
       include.aic = F,
       include.bic = F,
       include.loglik = F,
       include.groups = F,
       stars = c(0.01, 0.05, 0.1),
       table = F,
       custom.gof.rows = list("Year FE" = c('NO', 'NO', rep('YES',9), 'NO', 'NO'),
                              "Rating FE" = c(rep('NO', 3), rep('YES',8), 'NO', 'NO'),
                              "Issuer FE" = c(rep('NO', 8), 'YES', 'YES', 'YES', 'NO', 'NO'),
                              "Issuer$\\times$Year FE" = c(rep('NO',11), 'YES', 'YES'),
                              "Pseudo $R^2$" = sapply(m, function(x) {1-x$deviance/x$null.deviance}))
)

# =============================================================================================
# TABLE 3 
# NOTE: this is marginally different from the published version due to a few erroneously 
# duplicated observations in the original workfile (only in the reporting of T3, there is no 
# mistake elsewhere). The difference is mostly limited to the first decimal digit and does not 
# affect any of the results. T3 produced by this replication file is the correct one.
# =============================================================================================

# Nobs by money-ness and call decision
tab <- as.data.frame(dcast(sub_event_df, 'event_call ~ itm_mkt_pre')[,-1])
tab[3,] <- tab[2,]/(tab[1,]+tab[2,])*100
tab[4,] <- tab[1,]+tab[2,]

tab2 <- as.data.frame(dcast(sub_event_df[gfc==F], 'event_call ~ itm_mkt_pre')[,-1])
tab2[3,] <- tab2[2,]/(tab2[1,]+tab2[2,])*100
tab2[4,] <- tab2[1,]+tab2[2,]

tab3 <- as.data.frame(dcast(sub_event_df[postgfc==T], 'event_call ~ itm_mkt_pre')[,-1])
tab3[3,] <- tab3[2,]/(tab3[1,]+tab3[2,])*100
tab3[4,] <- tab3[1,]+tab3[2,]

tab <- cbind(tab, tab3)
remove(tab2, tab3)
tab <- rbind(tab[4,], tab[1:3,])

rownames(tab) <- c('No. bond-events', '~~no. no-calls', '~~no. calls', '~~\\% called')

addtorow = list()
addtorow$pos = list(-1, nrow(tab))
addtorow$command = c("\\hline \\hline  & \\multicolumn{3}{c|}{Full sample} & \\multicolumn{3}{c}{Post-GFC} \\\\\n",
                     "\\hline")

digimat <- matrix(0, nrow = nrow(tab), ncol = ncol(tab)+1)
digimat[4,] <- 1

print(xtable(tab, auto = F,
             align = 'l|ccc|ccc',
             digits = digimat),
      file = 'T3.tex', floating = F,
      include.colnames = T,
      include.rownames = T, size = 'footnotesize', sanitize.text.function=function(x){x},
      hline.after = c(0,nrow(tab)), add.to.row = addtorow)

# =============================================================================================
# (5) Return plots
# =============================================================================================

# dynamic per indexed period
toplot_df <- function(var, dataset, what, transform, freq) {
  
  if (what == 'mean') {
    
    df <- dcast(dataset[, mean(get(var), na.rm = T), .(event_id, get(freq))],
                'event_id ~ get', value.var = 'V1')
    
  } else if (what == 'weight.mean') {
    
    df <- dcast(dataset[, weighted.mean(x = get(var), w = vlm, na.rm = T), .(event_id,  get(freq))],
                'event_id ~ get', value.var = 'V1')
    
  } else if (what == 'last') {
    
    df <- dcast(dataset[, last(get(var)[!is.na(get(var))]), .(event_id,  get(freq))],
                'event_id ~ get', value.var = 'V1')
    
  }
  
  if (freq == 'week') {
    names(df)[-1] <- gsub('-', 'neg', paste0('w', names(df)[-1]))
    setcolorder(df, c('event_id', paste0('wneg', 9:1), 'w0'))
  } else {
    setcolorder(df, c('event_id', 'prenotice', 'notice', 'postnotice'))  
  }
  
  colnames(df)[2] <- 'base'
  
  if (transform == 'geom') {
    
    df <- df[, lapply(.SD, function(x) {100*(x/base-1)}), event_id]
    
  } else if (transform == 'diff') {
    
    df <- df[, lapply(.SD, function(x) {100*(x-base)}), event_id]
    
  }
  
  df <- as.data.table(t(df[, lapply(.SD, function(x) {
    
    dat <- winsorize(x, 0.01, 0.99)
    if (sum(!is.na(dat)) < 2) {
      c(NA_real_, NA_real_)
    } else {
      summary(lm(dat~1))$coefficients[1:2]
    }
    
  }), .SDcols = names(df)[-1]]))
  
  colnames(df) <- c('mean', 'se')
  df[, `:=` (min = mean - 1.96*se, max = mean + 1.96*se)]
  df = cbind(1:nrow(df), df)
  names(df)[1] = 'period'
  return(df)
  
}

# plot two prices overimposed
dat <- toplot_df(var = 'vwap', dataset = eventdf_plot[itm_pre == T & event_call == F & year >= 2010],
                 what = 'mean', transform = 'geom', freq = 'week')
dat[, called := F]
dat2 <- toplot_df(var = 'vwap', dataset = eventdf_plot[itm_pre == T & event_call == T  & year >= 2010],
                  what = 'mean', transform = 'geom', freq = 'week')
dat2[, called := T]
dat <- rbind(dat[1:9,], dat2[1:9,]); remove(dat2)

p <- dat %>%
    ggplot(mapping = aes(x=period, y=mean, group=called, colour=called)) + 
    geom_ribbon(aes(ymin=min, ymax=max, fill=called), alpha=0.4, colour = NA) +
    scale_color_brewer(palette = 'Set1', guide = 'none') +
    scale_fill_brewer(palette = 'Set1', labels=c('Not called', 'Called')) +
    geom_line(lwd = 2) +
    geom_point(size = 3, shape = 21, stroke = 3, fill = "white") +
    theme_bw(base_size = 32) +
    theme(panel.grid.major = element_line(size = 0.5, linetype = 'dashed', colour = "lightgray"), 
          panel.grid.minor = element_blank(),
          plot.title = element_blank(), legend.background = element_blank(), legend.position = c(0.1, 0.9),
          legend.title = element_blank(),
          axis.text=element_text(size = 30), axis.title=element_text(size = 30)) + 
    geom_hline(yintercept = 0, linetype = "dashed", color = "black", lwd = 1.5) +
    geom_vline(xintercept = 5, linetype = "dotted", color = "gray", lwd = 1.5) +
    scale_y_continuous(labels = scales::number_format(accuracy = 0.1, decimal.mark = '.')) +
    labs(x = "Week", y = 'Cumulative price return to week -9, %') +
    scale_x_continuous(labels = seq(-9, -1, 1), breaks = seq(1, 9, 1), expand = c(0.05, 0)) +
    annotate("text", x = 5, y = 0.7, label = "Notice deadline", hjust = 1.1, size = 9, fontface = 'italic')

ggsave(plot = p, filename = 'F2.pdf', width = 16, height = 9)